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THE DETERMINATION OF A STABILITY INDICATIVE FUNCTION 


FOR LINEAR SYSTEMS WITH MULTIPLE DELAYS* 

By John D. Shaughnessy and Yasutada Kashiwagi 

SUMMARY 

A theoretical study is made of the stability of a class of linear differential- 
difference equations with multiple delays. A direct method for determining the exact 
stability boundaries for homogeneous, linear differential-difference equations with con- 
stant coefficients and constant delays is formulated. This formulation results in a sta- 
bility indicative function, depending on a single parameter, which determines the number 
of roots of the transcendental characteristic equation that have positive real parts. It is 
proved that the system is stable if and only if this function has a value of zero. 

A second-order system with delays in the velocity and position feedback terms is 
considered as an example, and the stability regions for this system are determined for a 
range of delays and coefficients. It is observed that introduction of a delay has a definite 
destabilizing effect on the system, and introduction of a second delay has a compounding 
effect to further reduce stability. However, this example clearly illustrates that certain 
combinations of delays can stabilize an unstable system. This phenomenon is discussed 
from a theoretical point of view. 


INTRODUCTION 

Time delays or retarded actions are present in modern problems of guidance and 
control and the dynamics of manned and unmanned space vehicles. The effect of such 
delays is the subject of this paper. These phenomena may occur in several different 
ways; for example, in remote control of distant space vehicles the communications delay 
can adversely affect the stability of the overall control system. (See ref. 1.) Time 
delays in the engine response of large jet transports can seriously affect the handling 
qualities of the aircraft. In manned systems, lagging commands caused by slow human 
response can cause a normally stable system to become unstable. Hypervelocity entry 
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vehicles can lose aerodynamic stability as a result of flow-field lags caused by the space- 
craft motions and ablation. (See ref. 2.) As a final example, combustion delays in rocket 
motors can lead to erratic or intermittent running, possibly terminating with an explosion 
of the motor. (See ref. 3.) References 4 and 5 are comprehensive bibliographies of 
works published before 1960 that pertain to time-lag systems. Reference 6 contains a 
general discussion of the existence and stability of solutions of differential-difference 
equations. 

As of this writing there is no straightforward or practical method for determining 
the stability characteristics of linear systems with multiple time delays. However, it 
was suggested in reference 7 that the stability indicative function for linear systems with 
a single delay, introduced in reference 8, could be extended to include systems with multi- 
ple delays. The subject of this report is the extension of the stability indicative function 
to include such systems. 

The mathematics of delay problems often results in the so-called differential- 
difference equations, which are similar in several respects to ordinary differential equa- 
tions; however, there are significant dissimilarities. It is easy to show that a linear 
differential-difference equation with constant coefficients is formally like an ordinary 
differential equation of infinite order. Thus, the characteristic polynomial for the 
differential-difference equation obtained by using the classical Laplace transformation 
will be of infinite degree and will necessarily have a countably infinite number of char- 
acteristic roots. Another important distinction between ordinary differential equations 
and differential-difference equations is in the initial conditions. For the differential- 
difference equations, initial functions rather than initial constants must be specified. 

As will be shown, the stability of a large class of linear time-delay problems may 
be obtained by determining whether any roots of the characteristic equation have positive 
real parts. Even though these systems have an infinite number of characteristic roots, 
it is possible to determine how many of these roots have positive real parts. In this paper 
a method will be formulated for finding the stability of these systems in terms of the time 
delays and the coefficients of the differential-difference equations. 

SYMBOLS 


A real, positive constant 

a n ,b n real constants 

fj,f 2 functions of o> defined by equations (49) and (50) 
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Hj,H 2 functions of s defined by equation (5) 
I stability indicative function 

j positive integer 

k positive integer 

L characteristic function 

Lr real part of L(cr+iw) 

Li imaginary part of L(cr+ici>) 

l positive integer 

M, m positive integers 

min minimum 

N, n positive integers 

p positive integer 

s Laplace transform variable 

Tp position time delay 

T v velocity time delay 

t nondimensional time 

x dependent variable 

z positive integer 

a particular value of a 

/3 particular value of co 
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real constant 


C 

0 set of critical values 0*’s 

Oq,Qi, • • . nondimensional time delays 

jjj 

6 q, 01 , • • • critical values of 0 q,0i, • • • 

a real part of s 

t time 

U real constant 

o> imaginary part of s 

co* critical value of co 

Mathematical notation: 

arg s argument of s 

£a,b] closed interval a = t = b 

£a,bj semiclosed interval a = t < b 

1 fi 

sgn x sign of x: x/|x| 

[xIp] set of x’s having property P 

x( n )(t) nth derivative of x with respect to t 
A Jacobian determinant 

e is an element of 

A bar over a symbol denotes the prescribed value. 
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STABILITY OF LINEAR SYSTEMS WITH MULTIPLE DELAYS 


The following form of homogeneous linear differential-difference equations with 
constant coefficients and constant delays will be considered: 


a N x^(t) + b N x^(t - + ajj^x^-^t) + b N _i x ^ N "^(t - 

+ a Q x(t) + bqx^t - 6 0 ) = 0 


+ • • 


Equation (1) may be expressed more compactly as 

N 

^ a n x^ n ^(t) + b n x( n )(t - 6^j 
n=0 L 


= 0 


( 1 ) 


( 2 ) 


The coefficients a n and b n are real constants and the constant delays 0 n are 
nonnegative finite values and are not necessarily distinct. By using the Laplace trans- 
formation, the characteristic equation of equation (2) is found to be 




= 0 


(3) 


where s = a + ico (cr and co are real variables). 

Pontryagin (ref. 9) has shown that if a^f = 0, b^ 0, and 0^ > 0, the system is 
unstable. In reference 10 it is shown that if a^ ^ 0, bjq 7^ 0, and > 0, the system 
may be unstable even though all characteristic roots of equation (3) have negative real 
parts. These two cases will not be considered; rather, an examination will be made of 
the case where ajq- ^ 0 and bj^ = 0, written as 


N N-l 

L(s) = ^ a n s n + ^ b n s n e"^ nS = 0 fa N £ 0) (4) 

n=0 n=0 


Because of the infinite-degree nature of equation (4), exact determination of the 
characteristic roots is obviously impossible; however, the number of roots with positive 
real parts can be determined. Theorem 1 of the appendix can then be employed to 
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determine the stability of the system. If stability boundaries are desired, the number of 
roots with positive real parts is determined as a function of the parameters of interest. 

Before entering the development of this idea, several existing stability criteria are 
considered. The Routh-Hurwitz stability criterion clearly cannot be used with 
differential-difference equations. The Nyquist criterion can be used with only limited 
success because, as pointed out in reference 8, there is a question regarding the order 
of the characteristic equation which may lead to contradictory results. Also, stability 
boundaries would be quite difficult to obtain. Pontryagin (ref. 9) gives necessary and 
sufficient conditions for stability of linear systems with delays, but the delays must be 
integral multiples of each other. Finally, the Lyapunov stability criterion is applicable 
to linear as well as nonlinear differential-difference equations (see refs. 6, 11, and 12); 
however, the resulting stability boundaries depend on the Lyapunov function used and 
oftentimes do not give the largest region of stability. 

It was pointed out above that, in general, each characteristic root of equation (4) 
cannot be located, but the number of roots having positive real parts can be determined. 
Further, by theorem 1 of the appendix, the system is stable if, and only if, this number 
is zero. Before developing the method for calculating the number of roots that lie in the 
right half of the s-plane, some important aspects of this modified root-locus technique, 
extended to delay systems and first presented in reference 8, should be mentioned. 

Considerable information can be realized by studying a periodic solution. When 
some of the roots (nonmultiple) of equation (4) lie on the imaginary axis and the rest lie 
in the left half-plane, a periodic solution will be obtained. These pure imaginary roots 
correspond to intersections of the root loci of equation (4) with the imaginary axis. The 
root loci are the plots of the variations of the roots of equation (4) with changes in a 
single parameter, say an arbitrary delay 0^, from zero to its maximum value, while all 
other delays are fixed. The values of the parameter that correspond to an intersection of 
a root locus with the imaginary axis will be denoted as critical values 0j£. Also a point 
(0,icu*) in the s-plane will be denoted as a critical point. It will be shown that, in general, 
the number of critical values will be zero or countably infinite and these critical values 
can be found by replacing s with ico* in equation (4). If there are an infinite number 
of critical values, say 1 < e k 2 < ' * • < 0 k l < * ' * ( w ^ ere the subscripts 1,2,. . ., 

I, . . . indicate the different critical values of 0^), then the stability of the system 
remains unchanged in the interval 0 = 0^ < 0£ ^ but is unknown for % = \ because, 

so far, nothing is known about the behavior of the root loci at the intersections at s = iw*. 

The statements regarding the behavior of the root loci near the critical points will 
now be considered in more detail. This analysis follows directly from reference 8. A 
system which is stable when 0^ = 0 is now considered, and it is supposed that there are 
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two critical values 1 < e k 2 for e k € [ 0 > A ] where A is a real, positive constant. 
Thus the system is stable for 0 = 6^ < 0^ but the stability is not known for 0^ = 0^ 

There are several ways in which the root locus might intersect the imaginary axis 

when 0 k = 1 and % = ^ 2 ^ or 0 k 6 [o>aJ. Figure 1 illustrates three possible 

cases. Just the upper half of the s-plane is considered, since the loci are symmetric 
with respect to the real axis. If the loci have the form in figure 1(a), the system is 
stable for 0 = 0^ = 0^ j and unstable for 0^ ^ < 0^ = A. In the case of figure 1(b), the 

system is stable for 0 = % = 0^ j and 0^ 2 = % = A but unstable for 0 k 1 < < 0 k 2' 

In figure 1(c), the system is stable for 0 = 0k = A. Even if all three cases have the same 
critical values 0j* ^ and 0j* these cases represent different stability situations for 

0k 1 < ^k = A- Thus it is clear that in order to investigate the stability of systems 

described by equation (4), it is not sufficient to find only the critical values of a parame- 
ter. If only the critical values 0^ 0^ • • •> l’ ‘ ‘ ' are known for the given 

system, then all that can be said is that the system has the same stability characteristics 
for 0 § 0k < 0^ i ^where 0j£ 1 < e k 2 < • ' • < ®k l < ' ' ’) as for ®k = Tbus it is 
entirely possible for delays to make stable systems unstable and unstable systems stable. 
In general, if the stability of a system with delays is to be determined for 0k e ^,°°J , 
further information such as root loci directions at the critical points must be obtained. 

Calculation of the critical values and determination of the directions of the root 
loci at critical points are now considered. The characteristic equation can be written in 
the form 


L(s) = Hi(s) - H 2 (s)e" e k s = o 


(5) 


where 


N k-1 N-l 

Hi(s) = ^ a n s n + ^ b n s n e 0nS + ^ b n s n e 0n ‘ 

n=0 n=0 n=k+l 


and 


H 2 (s) = -b^ 
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If simple pure imaginary characteristic roots exist, there is at least one real positive 
constant co* such that 


-i<£o>* 

L(io>*) = Hi(ico*) - H2(ia>*)e K = 0 


( 6 ) 


and 


dL(s) 

ds 


s=io>* 


^0 


where the superscript * indicates critical values as before. Equation (6) may be 
divided by H 2 (ico*). If H 2 (iw*) = 0, then Hq(io>*) = 0 and the stability of the system 
does not depend on 0^; this case is not of interest here. Thus 


H l( iCtJ *) = e -ieJ w * 

H 2 (io>*) 6 


(7) 


From equation (7), the critical values of 0 k satisfy 


Hj(iw*) 

H 2 (iw*) 


0 k = — arg 
co* 


Hi(iw*) 

H 2 (iw*) 


(8) 


In general, the function Hj(iw*) will contain trigonometric functions of o>*, and 
equation (8) will not necessarily yield an explicit expression for u>*. If an explicit 
expression for co* cannot be obtained, an iteration technique will locate the co* values. 
In general the maximum value of co* for the iteration will have to be determined sep- 
arately. The number of critical values of 0k is zero or countably infinite, as can be 
seen in equation (8), since the various possible values of 0 k w* differ by integral multi- 
ples of 2 it. If the number of critical values is zero, the delay 0^ does not affect the 
stability of the system, since the absence of critical values means that no roots lie on the 
imaginary axis for that range of 0 k . If there are an infinite number of intersections, the 
way in which the loci intersect the imaginary axis for different critical values of 0 k 
must be found. 

There are three possible ways a root locus can intersect the imaginary axis. These 
cases are illustrated in figure 2. In case (a) the loci cut from the right half-plane to the 
left half-plane. In case (b) the loci cut from left to right. In case (c) the loci do not cut 
across the imaginary axis but do touch it. 
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The problem of finding how the root loci intersect the imaginary axis is now con- 
sidered. Since 6y. is taken by itself, the results of reference 8 may be extended to this 
analysis as follows. The behavior of the root loci at critical points (i.e., how they inter- 
sect the imaginary axis) can be determined by finding both the sign of the derivative 


.Ml 


'■a 


de 


M k 


^0 


0 k= 0 k 


and the value of M^, where 


Mjj is defined as 


Mjj = mini m 


»(e£) t o 


( 9 ) 


or a 


.(m) 


in 


= 0 for m = 1, 2, . . M^. - 1 and cr^k) j ^ 0. (The superscript i 

parentheses denotes the mth derivative.) If a root locus cuts across the imaginary axis 
from left to right, then 


and 


sgn = +1 

Mjj = Odd number^ 




(10a) 


If a root locus cuts across the imaginary axis from right to left, then 


sgn 


c (Mk) (4) = - 


and 




(10b) 


Mjj = Odd number^ 

Finally, if a locus intersects but does not cut across the imaginary axis, 

Mjj = Even number 


(10c) 


The quantity is determined as follows. Substitution of s = a + ia> into 

equation (5) yields 


L(c;+ico) = Hi(o+ico) - H 2 (cH-iw)e ^k( a+ * w ) 
= L R (a,co,0k) + iLi(a,w,%) 


= 0 


( 11 ) 
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Then it is clear that the root loci satisfy 


L R (a,co,6ik) = 0 


(12) 


Ll(cr,w,0k) = 0 


(13) 


In order to find the derivative of a with respect to 0 k , equations (12) and (13) are 
differentiated with respect to % to yield 


dL R 

d0 k 


9Lr 9Lr da 0 Lr dco 

+ * — -I — — - 

90 k 9a d0 k dco d0 k 




dLj _ 9Li 9Li^a_ 8L I dco _ 
d0 k 90 k 9 a d0 k dco d0 k 


In terms of the Jacobian determinant, 


da _ 1 8 (l-R> L l) 

de k a 9(e k ,o)) 


and 


where 


dco 

d0 k 


A = 


8 ( L R> L l) 

9(a,w) 


1 8 ( L R> L l) 


A 9 (a, 

Ok) 

9 ( L R’ L 1 

i 

9(a,u>) 


8L R 

9 Lr 

9a 

dco 

9 Li 

9L! 

9a 

9 CO 


> 0 


(14) 


(15) 


(16) 


(17) 
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and 


8 ( l r> l i) e 

8 ( 0 k^) 



0l R 

0 Lr 

9^k 

da> 

9Lj 

9Lj 

8% 

9w 

9LR 

0 LR 

9g 

90 k 

9Lj 

9Lj 

To 7 

90k 


( 18 ) 


(19) 


In reference 8 it is shown that the inequality in equation (17) will hold as long as the prob- 
lem is well defined; that is, if a^ is not equal to zero and b n 0 n (n = 0, 1, . . N - 1) 
are not all equal to zero. Thus, except for the case of branch points (referred to as crit- 


ical points in complex variable theory) where 


dL(s) 

ds 


j 

= 0, the derivatives and 

s=ico* 


— — given by equations (15) and (16) will exist. 

de k 


If i*£. = 0 at 0 k = 0 k ( note vertical tangent to root loci at critical points in 
de k 

figs. 2(a) and 2(b)), then M k > 1 and higher derivatives must be found in order to 
describe the root locus behavior at critical points. These higher derivatives are found 
from the following expressions: 


d m q 




+ 


_9_ 
9 CT 




8 

^ D 

T*-l 

1 

s 

\ dw 

9w 


|d0 k 


( 20 ) 


d m w 



9 

90 k 



+ 



dg 

d0 k 


+ 



dco 

d0 k 


( 21 ) 


for m > 1. 
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which are neces- 


Thus the value of and the sign of the derivative (^6^j 

sary for determining the behavior of the root loci at critical points may now be found by 
using equations (9), (15), and (20). The critical values of % are found by using equa- 
tions (8), as previously discussed. 


STABILITY INDICATIVE FUNCTION FOR LINEAR SYSTEMS 
WITH MULTIPLE DELAYS 


By using the results of the preceding analysis, the stability indicative function can 
now be formulated for linear systems with multiple delays whose characteristic equation 
has the form of equation (4), repeated here for convenience: 

( a N ^ °) 


N N-l 

L(s) = a n s n + ^ b n s n e" e n s = o 


n=0 


n=0 


It is noted that theorem 2 of the appendix insures that the stability characteristics 
of the system L(s) J ^_q = 0 are, in the limiting case, those of L(s) = 0 for positive 


values of 0^ « 1. If values of 0 n (where n = 0, . . ., N - 1, and n ^ k) are fixed and 
0^ is a parameter, and if the number of roots in the right half-plane for equation (4) with 
0^ = 0 is known, then the number of roots in the right half-plane for equation (4) can be 
determined as a function of 0^. This follows from theorem 2 together with the fact that 
the root loci are continuous. The remaining problem is that, in general, the number of 


roots of 


L(s) 



that lie in the right half-plane is not known for the case of multi- 


ple delays. A technique is now developed for determining this number. 


Initially all delays are set equal to zero, and the number of roots of the resulting 
Nth-degree polynomial having positive real parts and nonnegative imaginary parts is 
readily found. At this point one parameter, say 0 q, is introduced and is allowed to vary 
from zero to its prescribed value 0g. If critical values of 0g are found (by using 

eqs. (8)) in the interval 0 = 0q = 0q, say % 1 < 2 < • * • < ^0’ then starting with the 

smallest value, 0 q p Mq ^ and sgn o-(^0,l)^* are determined as previously dis- 
cussed. Then, by using equations (10), it is determined whether a locus of characteristic 
roots of this system goes into or comes out of the left half-plane as 0q is increased 
beyond 0 q j. Because of symmetry, only the upper half-plane need be considered. The 
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next critical value, 0g 2 , is considered, and Mg 2 and sgn cr °’^0g 2J are found as 
before. Again equations (10) are used to determine the behavior of the root locus at 

j|c 

0q = 0q 2 . After this procedure, it will be known precisely how many roots lie in the 
upper right half-plane when 0q = 6q. The next parameter, say Q\, is considered with 
6*0 = 0q and all other delays equal to zero. Now 6\ is allowed to vary from zero to 
01 and each critical value is taken account of as before, so that when = 0^, the num- 
ber of roots in the upper right half-plane will be known. This procedure is continued 
until only 0 ^ remains. 


At this point, the number of roots of L(s) 


0 k -° 


= 0 that have positive real parts 


and nonnegative imaginary parts (i.e., that lie in the upper right half-plane) will be known. 
The stability characteristics of equation (4) can now be expressed in terms of To do 

this, is allowed to vary from zero to its maximum value, taking account of each crit- 


ical value and finding and sgn 


* (l %) 


as above. If there are no roots of 


L(s) 


0 k =O 


0 in the first quadrant of the s-plane, then by theorem 1 the system is stable 


when 0k = Tli c system remains stable as 0k increases from zero until the first 


critical value is reached where M k is odd and is positive; beyond this point 

the system becomes unstable. If more roots pass into the first quadrant as 0k is 
increased further, an account of the number is kept. The system remains unstable until 
all these roots leave the first quadrant as 0k becomes larger still. On the other hand, 


if some of the roots of L(s) 


0 k =O 


= 0 do lie in the first quadrant, then by theorem 1 the 


system is unstable for 0k = 0. As before, the system remains unstable so long as there 


are any roots having positive real parts and nonnegative imaginary parts. Thus by keeping 
account of the number of roots in the first quadrant as 0k is increased from zero, the 
stability characteristics of the system are determined in terms of the parameter 0k. 


This concept is now formulated into the "stability indicative function. " As indicated 
above, all delays are set equal to zero in equation (4) to form the Nth-degree polynomial 


L(s) = L(s)| 


%-l _6l N-2-- • •-01-%-° 


= 0 


N-l 

= a N s N + ^ (a n + b n )s n = 0 
n=0 


( 22 ) 
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Then z, the number of roots of L(s) = 0 that lie in the first quadrant of the s-plane, 
may be found by classical techniques. Without loss of generality it is assumed that 
0 q / 0. Then 0-±, . . ., 0^-1 are set e( l ual to zero to form 


L 0 (s) = L(s) 


0 N-l _0 N-2- 


6 ^— 9 0 


= 0 


( 23 ) 


Next Zq, the number of roots of Lq(s) = 0 that lie in the first quadrant when 9q = 9q, 
is found. Since at each step in this procedure the number of roots in the first quadrant 
will be known initially, and since the loci are continuous as mentioned above, the "stability 
indicative function" for systems with a single delay introduced in reference 8 may be used 
to find zq, zj, and so forth. 

Now zq is equal to z plus any roots of Lg(s) = 0 that pass into the first quad- 
rant of the s-plane minus any that leave as 0q goes from zero to 9n. As explained 

(Mg) { *\ 

earlier, the value of Mq and the sign of cr v '(0 q) together determine the behavior of 
the root loci at the critical points. Thus, the stability indicative function for 6q is 
given by 


where 


I(0 O ) = z + 


0Q € ©0 


so that 


0 Q = 


o ^ e; s 0 O 


(24) 


(25) 


Zq = l(e 0 ) (26) 

In order to clarify further the meaning of I^0gj and I^0 q) , the following example 
is considered. It is assumed that there exists a system having z = 1 and the first five 
critical values of 9q such that 0 < 0q j < 6q 2 < 3 < ^0 4 < 5 < ^0* ^ b e h av i° r 

of the root loci is described by figure 3(a), then the function l(0g) would have the form 
shown in figure 3(b). It is clear that 1(0) = z = 1 and I^0 q^ = Zq = 0. Thus it is seen 
that the delay 9 g has made this normally unstable system stable. In this example the 
root loci are shown only for clarity; it must be remembered that in applying the stability 
indicative function no attempt need be made to plot the root loci. Instead, the stability 
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indicative function is the only quantity of interest and may be plotted or tabulated, for 
example. 

The stability indicative function for systems with multiple delays can now be 
developed. Another nonzero delay, say 0j with prescribed value 0j, is considered. 
Here Lq(s) = 0 is formed as 


Li(s) = L(s) 


0 N-l~ 0 N-2 = 


•~ 0 3 -0 2 -0, e O -0 O 


(27) 


and Zj, the number of roots of Lq(s) = 0 that lie in the first quadrant of the s-plane, 
is found. The value of zj is determined from the stability indicative function for 
given by 


I(6>i) = z 0 + 



(28) 


where 





so that 


Z 1 



(29) 


(30) 


The next delay, say 0 2 , is considered. Here L^s) = 0 is formed with 0 q = 6 q 

and 0j = 9 \) and Z 2 , the number of roots of L< 2 (s) = 0 that lie in the first quadrant of 
the s-plane, is found by using the stability indicative function. This procedure is con- 
tinued, considering one delay at a time, until 0k, the given parameter, is the only one 
remaining. 


At this point the number of roots of L(s) 


0 k =O 


= 0 that have positive real parts and 


nonnegative imaginary parts is known. Then z is defined as the number of character- 


istic roots of L(s) 


0 k =O 


0 whose real parts are positive and imaginary parts are non- 


negative. Now z can be expressed in terms of the above results as 
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Stability criterion .- If the stability indicative function for the parameter 0^ is 
constructed as 


If^k) = z + X, 


M k 


sgn a 


(Mk) 


k) 


e k e ©k 


(32) 


then the system whose transcendental characteristic equation is L(s) = 0 (eq, (4)) is 
stable if 


l(0k) = 0 (33) 

and unstable if 

I(0k) > 0 (34) 

Proof: By using equations (10) together with theorem 2 of the appendix, it is shown 
that the function 1(0^) given in equation (32) indicates the number of characteristic roots 
of equation (4) whose real parts are positive and imaginary parts are nonnegative for a 
given value of 0jj. Hence, from theorem 1 of the appendix, if l(0k) = 0 the system is 
stable, and if 1(0^) > 0 the system is unstable. 

STABILITY OF A SECOND-ORDER LINEAR SYSTEM WITH MULTIPLE DELAYS 

In order to illustrate the application of the stability indicative function for linear 
systems with multiple delays, the following differential-difference equation is considered: 

d - 2 - x ^ + 2 ?o dx . ( - 7 ^: .. Ty). + n 2 x(V - t p ) = o (35) 

The position delay T p is assumed to have a prescribed value and the stability 
boundaries of this system are to be found in terms of the velocity delay T v and the 
parameter To be compatible with the results of reference 8, which includes the 
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stability boundaries of equation (35) for the three cases where T v = T p , T v = 0, and 
Tp = 0, equation (35) is normalized as follows: 

Let 


t = Or 

Then equation (35) becomes 


0 O = T p £2 


— Tyf2 




dt 


dt 


( 36 ) 


(37) 


The corresponding characteristic equation is immediately obtained as 

L(s) = s 2 + 2£se~^l s + e s _ q (38) 


According to the theory, the first step in finding the stability characteristics of this 
system is to set the delays 6g and equal to zero and find the roots of the resulting 
second-degree polynomial. With 8q = = 0 there results 


L(s) = s 2 + 2£s + 1 = 0 


whose characteristic roots are 


s = 



(39) 


(40) 


It is clear that for £ > 0 there is no root in the first quadrant of the s-plane. Thus, 
z = 0 for £ = 0. 

Since &i is the normalized parameter of interest, the next step is to set 
equal to 0 and obtain 


Lg(s) = s 2 + 2£s + e 


-dnS 

0 = 0 


(41) 


By using equations (8) the critical values of 6 q and co for equation (41) are found to be 

(42) 


co* =\-2? 2 + y 4? 4 + 


g* _ 1 

y o — * 

u CO* 


tan 


-1 2 £' 


+ 2p7T 


CO 


(p = 0, 1, 2, . . .) (43) 
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The value of Mg and the sign of the derivative cr (6^) are now found by using equa- 
tions (9), (15), (16), and (20). From equations (12), (13), and (41), 


“ 0 0" 

L^(a,o),0o) = a 2 - a)2 + 2£a + e 0 cos QqU 




-0 n cr 


Lj(a,aj,0g) = 2co(a + £) - e 0 s i n gg W 


(44) 


J 


Then 



d0g 



1 8 ( L R» L l) 

A 8(0g,a>) 


e 0 =e% 



-2cj* - 2 £ a >* 0 * 
2£ - e*(w*) 2 


(co*) 

A 


4£ 2 + 2(cu*) 


2 


(45) 


Clearly 


da 

d0 o 


> 0 


8 <r e o 


so that Mg = 1 and sgn a W(4) = +1. Then for £ > 0 it is found from equations (24) 
and (26) that 


z 0 = 0 


for 





and 


z o = P 
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for 


r 


or 


U-1*L) + 

, cd*/ 


2(p - 1 )tt 


*o = 


Cl) 11 


tan -1 — M + 2pw 

CD*' 


Now the parameter 0^ is considered, with 8q = 6q = TpO. 

Lj(s) = s 2 + 2£se 6 1 S + e -6 0 S = 0 


By using equations (8), the expressions for the critical values of 
as 


(P=l, 2,. . .) 

Equation (38) becomes 
(46) 

6 1 and cd are found 


and 


(cd*) 4 - 4£ 2 (cd*) 2 - 2(a>*) 2 cos 0qcd* + 1 = 0 


(47) 



tan 


_1 (cd*) - COS 0qCD* 

sin "0QCD* 



(p = 0, 1, . . .) (48) 


The critical values of cd* are found for different values of £ by iterating equation (47) 
from cd = CD min to cd = cD max . Approximate values for CD min and cD max are found 
from equation (47) by solving for cos 0qcd* and finding the intersections of 


fj(CD) = cos 8qco 


with 


f 2 (cD) = 


cd 4 - 4^2 w 2 + i 
2cd 2 


(49) 


(50) 


by a graphical technique. 


The value of and the sign of the derivative 

From equation (46), 


ive oM 


0) 


are now determined. 


Lr(<7,cd,0i) = ct 2 - cd 2 + 2 %oe ®l a cos 0jcd + 2£cDe 4 sin 8\ cd + e ^°cos 0qCd (51) 
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and 


Li(ct,o), 01 ) = 2 aw - 2£ae ^l^sin 6\u> + 2£a>e -0 l°cos - e~^O a sin Oqw 


(52) 


Then 



(a>*) 4 + (cu*)^ 0 q sin 0qO>* - 1 


(53) 


Equation (53) is evaluated by using the given value of 0q and the critical values of w* 
found from equations (44). 

The stability boundaries of the system described by equation (35) may now be found 
by setting the stability indicative function for 0 ^ equal to zero, where 


1(01) = Zq 



sgn 


da 

d6>i 




(54) 


j _ 

Here it is assumed that M = 1; if = 0, then M / 1 and higher derivatives 

d$i 

h= e i 

must be found as discussed earlier. 

Equations (47), (48), (49), (50), and (54) were solved for a range of values of 0g 
and in this example, for 9 q between 0 and 2.6 and for £ between 0 and 1.2. The 
results of these computations are presented in figures 4 to 8 in the ($i,£)-plane with 
as a parameter. 

Figure 4 gives the stability regions of equation (37) for 0 q = 0. For small values 
of the system is stable. For small values of £ the stability of the system alter- 
nates as 0i is increased; however, at some value of 0i the system becomes unstable 
and it remains unstable for larger values of 0 1 > Note that the stability regions continue 
to repeat but diminish in size as 0i increases without limit. These results correspond 
to those obtained in reference 8. 
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Stable 


Unstable 





Figure 5 shows the stability boundaries of equation (37) for 0 q = 0.05. These 
boundaries are similar in some respects to those obtained when 0q = 0, but for values 
of £ less than approximately 0.025, the system is unstable for all values of 6\. It is 
seen that each region of stability for this case is located within each region of the pre- 
vious case. However, the stability regions vanish for sufficiently large values of 6\. 

The results of this example are slightly more useful when presented in the form of 
a topological plot of the stability regions. This is done for the three regions of figure 4 
in figures 6, 7, and 8. The coordinate scale factors are modified to allow enlarged plots. 
Figure 6 is probably of most practical value since the time delays will usually be small 
and will fall in the range of the variables considered. The area between each curve and 
the abscissa represents a region of stability for the particular value of 6q shown. It 
can be seen that when 6\ - 1 and £ = 1 the system is unstable for 0 = 0 q = 1.8 but 
becomes stable for 2.0 = ~6q = 2.4; that is, introduction of the second delay 0g has sta- 
bilized the normally unstable system. Figure 7 shows the destabilizing effects of the 
second delay 6q and indicates that the regions of stability vanish completely for 0q 
greater than approximately 0.3. Figure 8 shows a further reduction in the maximum 
value of 6q for which the system can be stable. Finally, it is noted that if the system 
is unstable for a given set of values in the (£,0i,7?o)-space, say (0.1, 4. 0,0.1), the unstable 
system can be made stable by increasing Oi to about 6; but increasing 0\ to the range 
shown in figure 8 or beyond will not stabilize the system. 

CONCLUDING REMARKS 

A theoretical study has been made of the stability of a class of linear differential- 
difference equations with arbitrary multiple delays. A direct method for determining the 
exact stability boundaries for homogeneous, linear differential-difference equations with 
constant coefficients and constant delays was formulated. The formulation resulted in a 
stability indicative function depending on a single parameter. This function determines 
the number of roots of the transcendental characteristic equation that have positive real 
parts for a given set of coefficients and delays, and it is proved that the system is stable 
if, and only if, this function has a value of zero. 

The stability indicative function for such systems compared favorably, in terms of 
range of application and ease of implementation, with other stability criteria. For exam- 
ple, Pontryagin’s criterion is limited to systems with integral multiple delays and the 
Lyapunov method generally does not give the largest region of stability without an exhaus- 
tive search for the proper Lyapunov function. A second-order system with delays in the 
velocity and position feedback terms was considered as an example, and the stability 
regions for this system were determined for a range of delays and coefficients. It was 
observed that introduction of a delay has a definite destabilizing effect on the system, and 
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Figure 6.- Stability boundaries of the system x(t) + 2£x(t - 0i) + xft - 0 q) = ° several values of 0 q w ' th 0 1 between 0 and 1.6. 









introduction of a second delay has a compounding effect to further reduce stability. How- 
ever, this example clearly illustrates that certain combinations of delays can stabilize an 
unstable system. 

Langley Research Center, 

National Aeronautics and Space Administration, 

Langley Station, Hampton, Va., September 20, 1968, 

124-07-02-1 6-23. 


30 



APPENDIX 


THEOREMS FOR THE DEVELOPMENT OF THE STABILITY INDICATIVE FUNCTION 

Theorem 1 

The homogeneous linear system whose characteristic equation is 

N N-l 

L(s) = a n s n + ^ b n s n e ^ nS = 0 
n=0 n=0 

(where the coefficients a n and b n are real constants and the constant delays 0 n are 
nonnegative finite values and are not necessarily distinct) is asymptotically stable if, and 
only if, all the characteristic roots have negative real parts. 

Proof.- See reference 11. 



Theorem 2 

For a sufficiently small positive value of 0k with bk ^ 0, the finite roots of 
L(s) = 0 can be made arbitrarily close to the finite roots of L(s)| g q = 0, and there 

exist an infinite number of roots whose real parts are negatively infinite. 

Proof .- Express the roots of L(s) =0 as s = a + ico, consider the following 
assumptions, and prove the theorem by contradiction. 

Case 1: lim o = a and lim to = /3 where a and 0 are finite constants. 

0k~*'O 0k"*O 

Thus, 


N 


k-1 N-l 

^ 1 r. -0„S ST 1. „ -0„S 


lim L(s) ='y * a n s n + bkS k + /> 'b n s n e nS + ^> 'b n s n e nS = 0 
? k“° n=0 n=0 n=k+l 


0 k -O 

Hence, we see that the roots of L(s) = 0 are close to the roots of L(s) 
small values of 0k- 


0 k =O 


= 0 for 
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Case 2: lim a = a and lim o> = °° where a is a finite constant. Then 

e k -0 0 k -o 


lim |e”^n s | _ ii m e ~ d n° ~ e ~ d n a 
e k-° 



lim ± = 0 

\Jcr 2 + (jj2 


Divide both sides of L(s) =0 by s N £ 0 to get 


N 

I 

n=0 


U N-l 

a n s n- N + ^ b n s n "N ( 
' n=0 


« e " 0 n s = 0 


N-l / \ 

a N + ^ ( a n + b n e~ e n s J s n-N _ 0 
n=0 

If case 2 holds, then each term in brackets is finite and a^ = 0. This result con- 
tradicts the definition of L(s) = 0; hence, case 2 does not hold. 

Case 3: lim a = +<*>. Then 

0k-O 


lim je“ 0 n s | = ii m e -0 n cr = o 
0 k -*-O 0 k —O 


and 


Iim - 
0 k -O s 


lim - = 0 

0 k~° \jo 2 + co 2 


Again, this result requires that a^ = 0. Thus case 3 does not hold, and the only 
other possibility wherein ajg ^ 0 is that 
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e k~° 
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